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Abstract 


We  consider  the  following  four  problems  for  a  set  S'  of  fc  points  on  a  plane,  equipped 
with  the  rectilinear  metric  and  containing  a  set  R  of  n  disjoint  rectangular  obstacles  (so 
that  distance  is  measured  by  a  shortest  rectilinear  path  avoiding  obstacles  in  R):  (a)  find 
a  closest  pair  of  points  in  S,  (b)  find  a  nearest  neighbor  for  each  point  in  S,  (c)  compute 
the  rectilinear  Voronoi  diagram  of  S,  and  (d)  compute  a  rectilinear  minimal  spanning  tree 
of  S.  We  describe  0((n  +  k)  log(n  +  k))  time  sequential  algorithms  for  (a)  and  (b)  based  on 
plane-sweep,  and  the  consideration  of  geometrically  special  types  of  shortest  paths,  so-called 
z-first  paths.  For  (c)  we  present  an  0((n  -t-  k)  log(n  -|-  k)  log  n)  time  sequential  algorithm 
that  implements  a  sophisticated  divide- and- conquer  scheme  with  an  added  extension  phase. 
In  the  extension  phase  of  this  scheme  we  introduce  novel  geometric  structures,  in  particular 
so-called  z-diagrams,  and  techniques  cissociated  with  the  Voronoi  diagram.  Problem  (d)  can 
be  reduced  to  (c)  and  solved  in  0((n  -I-  fc)log(n  -1-  ^;)logn)  time  as  well.  All  our  algorithms 
are  near-opiimalj  cls  well  as  easy  to  implement. 

Keywords:  Computational  geometry,  rectilinear  metric,  obstacles,  nearest  neighbors,  Voronoi  diagram,  min¬ 
imal  spanning  tree. 
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A  fundamental  problem  in  computational  geometry  is,  given  a  geometric  space  G  equipped  with 
some  metric  d,  that  of  computing  shortest  paths  in  G.  This  leads  to  several  proximity  problems 
when  one  is,  in  addition,  given  a  finite  subset  S  of  G. 

We  consider  the  case  when  G  is  the  Cartesian  plane,  together  with  an  obstacle  set  R  consisting 
of  n  disjoint  isothetic  rectangles  (i.e.,  with  sides  parallel  to  the  coordinate  axes),  and  the  metric 
d  (often  called  the  rectilinear  oi  Manhattan  metric)  is  defined  such  that,  if  p,  5  €  G,  then  d{p^q) 
is  the  Euclidean  length  of  a  shortest  rectilinear  path  (i.e.,  consisting  of  axes-parallel  segments) 
joining  p  and  q  that  does  not  intersect  the  interior  of  any  of  the  rectangular  obstacles  in  R.  Such 
geometric  spaces  arises  naturally  in  applications  such  as  VLSI  chip  design,  plant  and  facility 

*An  extended  abstract  appeared  in  Proc.  13th  Conf.  on  the  Foundations  of  Software  Tech,  and  Th.  Comp.  Sc., 
Bombay,  1993,  Springer- Verlag,  218-227. 
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^Supported  in  part  by  the  National  Science  Foundation  under  grants  CCR-9004346  and  IRI-9307506,  the  Office 
of  Naval  Research  under  grant  N00014-94-1-0284,  and  an  endowed  chair  supported  by  Hitachi  Ltd.  at  Faculty  of 
Engineering  Science,  Osaka  University. 
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layout,  robot  motion  planning,  and  urban  transportation.  These  spaces,  as  well  as  various 
generalizations,  have  been  investigated  extensively  for  efficient  sequential  [5,  6, 11, 13,  14, 18, 19] 
and  parallel  [2,  8]  algorithms  to  compute  shortest  paths. 

The  problems  considered  in  this  paper  are,  given  a  plane  G,  containing  an  obstacle  set  R  of 
n  rectangles  and  equipped  with  a  metric  d,  as  described  above,  and,  further,  a  set  S  of  k  points 
in  G  -  Ur^R  r,  to 

(a)  find  a  closest  pair  of  points  in  5, 

(b)  find  a  nearest  neighbor  for  each  point  in  5, 

(c)  compute  the  Voronoi  diagram  of  5,  and 

(d)  compute  a  minimal  spanning  tree  of  5. 

We  give  near-optimal  sequential  algorithms  for  all  four  problems  (in  fact,  each  algorithm 
may  be  sub-optimal  by  at  most  a  logarithmic  factor),  thus  resolving  questions  open  since  at 
least  1985  when  de  Rezende,  Lee,  and  Wu  [6]  gave  optimal  algorithms  to  find  shortest  paths  in 
such  a  space. 

Section  2  introduces  some  terminology  and  preliminary  results,  as  weU  as  the  notion  of  so- 
caUed  z-first  paths  (where  z  is  one  of  the  four  directions,  and  iy),  which  are  shortest  paths 
of  a  special  geometric  type. 

In  Section  3,  we  describe  an  0((n  +  k)  log(n  -f*  k))  time  algorithm  that  finds  a  closest  pair  in 
S  after  sweeping  the  plane  in  the  -t-cc  and  -\-y  directions  to  determine  shortest  x-first  and  y-first 
paths. 

Section  4  describes  a  similar  algorithm  with  the  same  time  bound  to  find  aU  nearest  neighbors 
in  5,  but  in  this  case  the  plane  is  swept  in  aU  four  directions,  ±x  and  ±y. 

In  Section  5,  we  describe  the  more  complicated  algorithm  to  compute  the  Voronoi  diagram  of 
S.  Our  algorithm  runs  in  0((n-l-  A;)log(n-t-  k)logn)  time,  which  is  quicker,  but  not  significantly 
so,  than  the  next  best  Voronoi  diagram  algorithm  that  we  are  aware  of  for  a  similar  geometric 
space:  the  algorithm  of  Mitchell  [14,  see  Theorem  2]  which  runs  in  0((n  +  A:)log^(n  +  k))  time. 
However,  what  may  be  of  more  interest  is  that,  while  Mitchell’s  algorithm  uses  a  “continuous 
Dijkstra”  method  of  propagating  a  “wavefront”'  from  each  point  of  5  as  a  source,  ours  is  quite 
different  and  based  on  a  divide-and-conquer  with  an  “extension  ph«ise”.  It  should  be  pointed 
out  though  that  Mitchell’s  method  allows  the  more  general  class  of  simple  polygons  as  obstacles, 
while  it  is  not  clear  if  our  methods  can  be  extended  beyond  the  class  of  rectangular  obstacles. 

Considering  another  geometric  space  with  obstacles,  Aronov  [1]  achieves  a  time  bound  iden¬ 
tical  to  ours  for  computing  the  Euclidean  Voronoi  diagram  of  k  points  in  an  n-sided  simple 
polygon.  His  overall  scheme  is  also  divide-and-conquer  with  an  extension  phase  but,  our  space 
being  dissimilar,  we  differ  significantly  in  the  implementation  of  the  scheme  and,  in  fact,  intro¬ 
duce  new  geometric  structures  and  methods.  In  particular,  we  exploit  the  special  geometry  of 
the  rectilinear  plane  to  implement  the  crucial  extension  phase  in  two  stages:  in  the  first  stage 
we  compute  “approximate”  Voronoi  extensions,  so-caUed  z-diagrams  (where  .z  is  one  of  the  four 
directions,  ±x  and  ±y),  and  then,  in  the  next  stage,  use  these  approximations  to  sweep  through 
while  tracing  out  the  boundaries  of  the  “exact”  Voronoi  cells. 

In  Section  6,  we  discuss  the  reduction  of  the  problem  of  computing  a  minimal  spanning  tree 
of  S  to  that  of  computing  the  Voronoi  diagram  of  S.  This  leads  to  an  0{{n-{- k)\og{n  + k)  log  n) 
time  algorithm  to  compute  a  minimal  spanning  tree,  implying  almost  linear  order  speed-up  over 
the  minimal  spanning  tree  algorithm  of  Wu,  Widmayer,  Schlag,  and  Wong  [18]  (which  runs  in 
0{klogk  -h  TO^logn)  time,  but  allows  the  more  general  class  of  rectilinear  convex  polygons  as 
obstacles). 

We  conclude  in  Section  7  with  a  discussion  of  the  near-optimality  of  all  our  algorithms  and 
related  open  questions. 

Throughout,  we  avoid  repeating  proofs  that  have  appeared  in  the  available  literature. 
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Figure  1:  The  shaded  rectangles  comprise  a  horizontal  blockage  between  s  and  g;  P  is  the  z-lirst 
path  from  p  to  q. 

2  Preliminaries 

We  shall  henceforth  assume  that  all  paths  are  rectilinear  and  avoid  intersecting  the  interiors  of 
rectangles  in  R.  A  path  in  G  from  a  point  p  to  a  point  q  that  is  monotone  in  the  x  direction  is 
called  an  x~path,  a  path  monotone  in  both  x  and  y  directions  is  called  an  xy-path,  and  similarly 
denote  paths  monotone  in  other  directions  (-a:,  -y,  etc.). 

A  horizontal  blockage,  if  one  exists,  between  a  pair  of  points  s  and  ^  is  a  subset  of  rectangles 
of  R  dispositioned  as  indicated  in  Figure  1  (in  which  case,  we  also  say  s  is  horizontally  blocked 
from  q).  More  precisely  (following  [6]): 

Assume  s  =  {s^iSy)  and  q  =  {qx^qy)  s-nd,  without  loss  of  generality,  that  Sx  <  qx 
Sy  <  qy.  For  a  rectangle  r  £  R,  let  minx{r')  and  maXx{r)  denote  the  smaller  and  larger  of  the 
x-coordinates  of  the  vertical  edges  of  r,  respectively.  Similarly,  define  miny{r)  and  TnaXy{r). 

Given  a  sequence  H  =  {ri,  ^2, . . . ,  r/j}  of  rectangles  of  R,  we  say  that  H  is  3.  horizontal 
blockage  from  s  to  g  if 

1.  Sy  <  minyiji), 

2.  maxy{r{)  <  mmy(r,+i)  Vi  =  l,2,...,/i-  1, 

3.  maxy{rh)  <  qy, 

4.  maXx(ri)  >  Sx  >  rninx{Ti), 

5.  maxx{ri)  >  rninx{TiJ^])  Vi  =  1,2, . .  .,h  -  1,  and 

6.  maxx{rh)  >  qx  >  minx{rh)- 

Similarly  define  vertical  blockage.  We  refer  to  [6]  for  a  proof  of  the  following: 

Proposition  1  ([6])  Between  two  points  p  and  q  there  cannot  be  both  a  horizontal  as  well  as 
a  vertical  blockage  (there  may  be  neither).  If  there  is  no  horizontal  blockage  and  p  is  left  of  q 
(i.e.,  the  x-coordinate  of  p  <  the  x-coordinate  of  q),  then  any  shortest  path  from  p  to  q  is  an 
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x-path.  If  there  is  no  horizontal  blockage  and  p  is  right  of  q,  then  any  shortest  path  from  p  to  q 
is  a  {—x)-path. 

Similar  results  hold  if  there  is  no  vertical  blockage.  ^ 

Shortest  paths  joining  two  points  are  never  unique  (except  in  the  trivial  case  when  the 
shortest  path  consists  of  only  one  segment),  but  in  the  following  proposition  we  introduce  and 
define  a  special  type  of  shortest  paths,  z-first  paths.,  with  a  geometric  property  that  almost 
always  makes  them  unique: 

Proposition  2  If  there  is  no  horizontal  blockage  between  p  and  q  and  p  is  left  of  q,  then  there  is 
at  least  one  shortest  x-path  from  p  to  q,  call  such  an  a: -first  path  from  p  to  q,  that  always  proceeds 
in  the  x  direction  unless  it  would  either  enter  inside  a  rectangle,  or  enter  a  region  between  each 
point  of  which  and  q  there  does  not  exist  an  x-path  (i.e.,  each  point  of  that  region  is  either 
horizontally  blocked  from  q  or  lies  to  the  right  of  q).  Thus,  the  x-first  path  from  p  to  q  makes 
a  turn  in  a  ±y  direction  only  when  it  either  hits  a  rectangle,  or  “risks”  loosing  x -monotonicity. 
See  Figure  I. 

Similar  results  hold  if  p  is  right  of  q,  or  if  there  is  no  vertical  blockage.  In  particular,  we 
have  similar  definitions  for  z-first  paths,  where  z  =  —x  or  ±y. 

Proof.  A  precise  iterative  procedure  for  drawing  an  x-first  path  from  p  to  g  is  as  follows: 

Say  the  source  p  =  (PxiPy)  destination  q  =  {qxiQy)-  Draw  the  path  from  p  in  the  x 
direction  until  (whichever  comes  first) 

1.  it  hits,  at  point  c,  the  left  edge  of  some  rectangle,  or 

2.  it  reaches  a  point  d  =  {a,Py)  such  that,  for  sufficiently  small  €  >  0,  {a-{-e,Py)  is  horizontally 
blocked  from  q,  or 

3.  it  reaches  the  point  e  =  {qx,Py)- 

In  case  1,  say  the  corners  of  the  edge  on  which  c  lies  are  a  and  b.  Draw  the  path  vertically 
from  c  to  a  if  d{c,  a)  -|-  d{a,  q)  <  d{c,  b)  -f  d{b,  q),  or  to  b  if  d(c,  a)  -f  d{a,  q)  >  d{c,  b)  -H  d{b,  q),  or, 
arbitrarily,  to  either  a  or  6  if  d(^c,a)  -f-  d{a,q)  —  d{c,b)  -f-  d(b,q)  (this  is  exactly  the  case  when 
the  x-first  path  is  not  unique).  See  Figure  1.  Repeat  the  drawing  procedure  with  the  current 
endpoint  (either  a  or  6)  of  the  path  as  the  new  source. 

In  case  2,  let  the  rectangle  r,  from  a  horizontal  blockage  between  (a  -f-  e,Py)  and  q,  be  the 
one  that  is  vertically  adjacent  to  (a  -t-  €,Py),  for  sufficiently  small  e  >  0.  Then,  clearly,  the  left 
edge  of  r  has  corners  a  and  b  with  x-coordinate  equal  to  a.  Draw  the  path  vertically  from  d 
to  the  more  distant  of  o  and  b.  See  Figure  1.  Repeat  the  drawing  procedure  with  the  current 
endpoint  (either  a  or  b)  of  the  path  as  the  new  source. 

In  case  3,  draw  the  path  vertically  from  e  to  q.  This  is  possible  as,  by  case  2,  we  never  reach 
a  point  that  is  horizontally  blocked  from  q.  See  Figure  1.  Exit. 

Clearly,  this  procedure  completes,  after  a  finite  number  of  turns,  an  x-path  P  from  p  to  q. 
It  remains  to  show  that  P,  which  we  caU  the  x-first  path  from  p  to  q,  is  indeed  shortest.  We 
shall  prove  this  by  induction  on  the  number  of  segments  of  P. 

Starting  the  induction  is  trivial.  Assume  inductively  then  that  P  has  n(>  1)  segments,  and 
that  all  x-first  paths  with  no  more  than  n  -  1  segments  are  shortest.  If  possible  let  Q  be  a  path 
from  p  to  g  that  is  shorter  than  P.  Then  P  and  Q  are  disposed  as  in  either  Figure  2(a)  or  2(b) 
(we  may  assume  without  loss  that  they  do  not  intersect). 

Assume  first  that  they  are  disposed  as  in  Figure  2(a).  Let  m  be  the  other  end  of  the  first 
segment  of  P  leaving  p  (this  segment  is  depicted  as  horizontal  in  Figure  2,  but  it  may  as  well  be 
vertical).  From  m  draw  the  xy-path  R  with  y-preferred  (i.e.,  the  path  which,  whenever  it  can 
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Figure  2:  Illustration  for  the  proof  of  Proposition  2. 


go  in  either  the  x  or  y  direction  without  intersecting  the  interior  of  a  rectangle,  chooses  the  y 
direction,  see  [6]).  As  m  is  not  horizontally  blocked  from  q,  R  intersects  Q  at  some  interior  point 
n.  Denoting,  for  example,  the  distance  along  P  from  ptoqhy  \P{p,  g)|,  we  have  by  assumption 

\Q{p,q)\  <  \P{p,q)\ 

^  \Q{p,n)\  +  \Q{n,q)\  <  \Pip,rn)\  +  \P{m,q)\  (1) 

As  P{p,m)  U  R[m,n)  is  an  xy-path,  it  is  a  shortest  path  from  p  to  ra,  so  that 

\P{p,Tn)\  +  \R(m,n)\  <  |i3(p,n)|  (2) 

From  (1)  and  (2)  we  have 

\P{p,m)\  +  \R{m,n)\  +  \Q{n,q)\  <  \Pip,m)\  + \P{m,q)\ 

\R{m,n)\  +  \Q{n,q)\  <  |F(m,g)|, 

contradicting  the  inductive  hypothesis  as  P{m,  q)  is  an  x-first  path  with  n  —  1  segments,  and 
proving  Q  cannot  exist.  If  P  and  Q  are  disposed  as  in  Figure  2(b),  an  exactly  similar  argument 
holds  after  drawing  R  from  m  as  the  a:t/-path  with  (— j/)-preferred.  O 

Comment.  The  notion  of,  for  example,  x-preferred  paths  in  [6]  is  different  from  ours  of  x-first 
paths  in  that  it  does  not  characterize  shortest  paths  between  given  pairs  of  points. 

The  following  separator  result  can  be  proved  following  [2]  (where,  in  fact,  an  n-processor 
logn-time  PRAM  algorithm  is  given)  with  straightforward  modifications,  also  see  Figure  3: 

Proposition  3  ([2])  In  0{n\ogn)  time  one  can  describe  an  x{-y)-path  P  (after  re-orienting 
the  coordinate  axes  if  necessary),  comprising  0{n)  segments  and  unbounded  in  both  directions 
with  the  first  and  last  segments  being  vertical  (P  is  imagined  to  start  from  a  point  at  infinity  in 
the  y  direction),  such  that  there  are  at  least  rectangles  of  R  on  either  side  of  P.  □ 


3  Finding  a  Closest  Pair 

It  is  worth  noting  at  the  outset  that  the  Bentley-Shamos  [3]  divide- and- conquer  scheme  cannot 
be  directly  applied,  as  it  is  no  longer  true  that  a  “circle”  of  radius  r  in  the  d-metric  around  a 
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Figure  3:  P  is  a  separator  as  in  Proposition  3. 


Figure  4:  Each  unlabeled  point  is  at  distance  less  than  r  from  p,  but  at  distance  greater  than  r 
from  each  other. 
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point  p  can  contain  at  most  0(1)  points  each  of  which  is  at  distance  at  least  r  from  the  other. 
See  Figure  4.  Instead,  we  exploit  properties  of  the  rectilinear  metric  for  a  very  different  method 
based  on  plane-sweep.  Our  idea  is,  while  sweeping  in  the  x  direction  for  instance,  to  gather 
information  about  shortest  x-first  paths.  It  should  be  pointed  out  that  the  algorithm  in  [6]  for 
the  single-source  shortest  path  problem  in  a  similar  setting  uses  plane-sweep  as  weU,  though,  for 
our  purposes,  we  need  to  maintain  considerably  more  information  and  data  structures  through 
the  sweep. 

Let  Cl  denote  the  set  of  corners  on  the  left  edges  of  rectangles  of  R.  Further,  as  a  simplifying 
device,  add  to  5  a  new  point  (-oo,0)  (we  shall,  in  fact,  assume  (-oo,0)  to  be  a  finite  point 
sufficiently  far  left  of  any  existing  member  of  5  U  Ci). 

Suppose  Lx  is  a  vertical  sweep-line  that  begins  to  scan  in  the  x  direction  starting  from  a 
position  just  right  of  (-oo,0).  Maintain  the  status  of  in  a  height-balanced  search  tree  (e.g., 
a  red-black  tree)  Tx,  such  that,  if  at  time  t  the  sweep-line  lies  on  Lx{t),  the  status  tree  at  that 
instant,  denoted  Tx{t),  represents  an  increasing  set  of  points  {-oo  =  ao,ai,. .  .,ai^  =  oo}  on 
LxH)  iLx{t)  is,  of  course,  imagined  to  be  a  copy  of  the  real  line  by  a  projection  of  the  y  axis), 
together  with  a  label  Bi  for  each  point  a,-,  1  <  i  <  ft,  such  that 

1.  For  each  i,  either  Bi  is  of  the  form  (p,  c),  where  p  G  5  and  c  £  S  U  Ci,  both  p  and  c  lying 
left  of  Lx{t)]  or,  Bi  is  of  the  form  r,  where  r  £  R. 

2.  If  Bi  =  (p,c)  then,  for  each  point  q  in  the  open  interval  (a^-i,  a,)  on  Lx{t),  p  is  the  closest 
to  q  of  the  points  of  S  from  which  there  is  an  x-path  to  q  (clearly,  all  such  points  must  lie 
left  of  Lx{i),  and  (-oo,  0)  is  one  such);  further,  there  is  a  unique  x-first  path,  say  P,  from 
p  to  q,  and  the  last  point  of  5  U  Q  at  which  P  turns  as  it  moves  from  p  to  g  is  c  (counting 
p  itself  as  the  first  and  possibly  only  such  turn  point).  Call  c  the  x-anchor  of  q  w.r.t.  p. 

3.  If  Bi  =  r  then  either  (aj_i,o,)  is  the  left  edge  of  r  (minus  endpoints)  or  lies  in  the  interior 
of  r. 

4.  The  sorted  set  {ao,  ui,  ■  •  • ,  uq}  is  minimal  with  respect  to  the  labeling,  i.e.  no  two  adjacent 
Oi  share  the  same  label. 

See  Figure  5  .  We  shall  often  refer  to  an  open  interval  as  being  associated  with  the  same  label 
as  its  right  endpoint.  Initially,  of  course,  Tx{0)  contains  only  the  points  — oo  and  oo,  with  the 
point  00  having  label  ((-oo,  0),  (-oo,  0)). 

In  addition  to  the  status  tree  Tx,  maintain  an  array  D,  indexed  by  SuCi,  such  that,  at  time 

•  if  g  is  left  of  Lx{t)  (i.e.,  if  Lx  has  already  swept  over  g),  the  entry  D[q]  contains  the  name 
of  the  point  which  is  closest  to  g  amongst  points  of  S,  different  from  g,  from  which  there 
is  an  x-path  to  g,  as  weU  as  the  distance  of  that  point  from  g,  and 

•  if  Lx{t)  is  strictly  left  of  g,  then  D[q]  =  oo,  denoting  ’no  information  available’  (as  a  special 
case,  mark  i?[(— oo,0)]  =  oo  throughout). 

As  Lx  moves  rightwards  the  following  three  types  of  events  occur  at  various  times: 

1.  Lx{t)  touches  a  point  of  S,  called  a  point  event. 

2.  Lx{t)  touches  the  left  side  of  a  rectangle,  called  a  left  event. 

3.  Lx(t)  touches  the  right  side  of  a  rectangle,  called  a  right  event. 
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Figure  5:  The  status  of  the  sweep  line  is  shown  at  five  positions:  (a)  just  before  the  point  event 
at  ps,  (b)  after  updating  for  the  point  event  at  ps,  (c)  after  updating  for  the  point  event  at  ps, 
(d)  after  updating  for  the  left  event  at  the  left  edge  of  tj,  and  (e)  after  updating  for  the  right 
event  at  the  right  edge  of  ry.  The  solid  circles  are  the  points  p,-  of  S,  the  hollow  circles  are  points 
of  the  status  tree,  the  n  are  the  rectangles,  and  the  c,-  are  corners  on  the  left  edges  of  rectangles. 
The  labels  are  shown  next  to  the  intervals,  and  label  b  =  ((— oo,0),C2). 


We  assume,  for  simplicity,  that  no  two  events  happen  at  the  same  time  (i.e.,  points  of  S  and 
vertical  edges  of  rectangles  of  i?  all  have  distinct  ^-coordinates).  Removing  this  restriction  is  a 
minor  technicality.  Sort  the  points  of  S  and  the  vertical  edges  of  rectangles  of  R  by  x-coordinate, 
in  0{{n  +  k)log{n  -|-  k))  time,  to  obtain  the  schedule  of  events.  The  following  is  not  hard  to 
prove: 

Lemma  1  Tx{t)  does  not  change  between  successive  events. 

Proof.  Follows  by  elementary  geometric  arguments,  given  the  rectilinearity  of  the  metric.  □ 

Next,  we  describe  the  procedures  to  update  D  and  Tx  at  each  type  of  event. 

Point  event:  At  a  point  event  at  time  Lx(t)  touches  a  point  of  5,  say  q.  Assume  the  previous 
event  took  place  at  time  t\ 

By  searching  Tx{t^)  locate  q  (or,  more  accurately,  its  y-coordinate)  in  the  array  ao, . . 

Say  ai-i  <  q  < 

If  the  label  Bi  =  (p,  c)  (it  cannot  be  an  r  G  R),  then  p  is  a  nearest  one  to  q  (we  have  a  choice 
when  g  =  a^)  amongst  points  of  5,  different  from  g,  from  which  there  is  an  x-path  to  g;  and  c 
is  the  x-anchor  of  g  w.r.t.  p.  Further,  the  distance  along  a  shortest  x-path  from  p  to  g  can  be 
determined  as  the  sum  of  the  distance  of  c  from  p  (which  can  be  read  from  array  D)^  the  length 
of  the  perpendicular  from  c  to  Lx{t),  and  the  distance  of  g  from  the  base  c  of  this  perpendicular 
(which  lies  in  (at_i,ai)).  In  fact,  the  function  di  :  [ai^i^ai]  3i,  the  reals,  giving  distances  of 
points  from  p  along  x-paths,  attains  a  minimum  at  c  and  linearly  increases  with  gradient  one 
on  either  side  of  that  base. 

Updating  D:  Update  the  entry  D[q]  with  the  name  of  p  and  its  distance  from  g. 

Updating  T^:  To  update  Tx{t^),  a  new  interval  with  label  (g,  g)  must  be  created  that  contains 
g  and  extends  on  either  side  of  g  to  contain  points  of  Lx{t)  that  are  now  closer  to  g  than  they 
are  to  points  of  S  strictly  left  of  Lx{t),  considering,  of  course,  only  distances  measured  along 
x-paths.  Determixiixig  the  extent  of  this  new  interval  is  straight-forward:  simply  proceed  in 
either  direction  from  g  along  Lx{t)  comparing  the  distance  from  g  to  the  distance  from  the 
currently  known  nearest  point  of  5  (as  given  by  the  distance  function  dj  of  the  interval  one  is 
currently  inside),  and  stopping  only  if  the  two  distances  become  equal  or  the  side  of  a  rectangle 
of  R  is  reached.  This  will,  of  course,  result  in  the  partial  or  full  deletion  of  some  vertices  of 
Tx{t').  Updating  Tx{t')  consists,  therefore,  of  inserting  the  two  endpoints  of  the  new  interval 
and  deleting  old  vertices  that  lie  inside  the  new  interval.  The  upper  endpoint  of  the  new  interval 
will  have  label  (g,g),  while  the  lower  endpoint  will  have  the  label  of  the  interval  of  Tx{t^)  in 
which  it  lies. 

See  Figure  5. 

Left  event:  At  a  left  event  at  time  t,  Lx{t)  touches  the  left  side,  say  [7'i,r2],  of  some  rectangle 
r  G  i2,  and  we  imagine  r  to  become  active  as  an  obstacle  intersecting  Lx{t). 

Updating  D:  Update  D[ri]  and  D[r2]  by  searching  Tx{t')  to  locate  vi  and  r2  in  qq,  • 
and  performing  appropriate  distance  computations. 

Updating  Tx'.  Updating  Tx{if)  consists  firstly  of  inserting  ri  and  r2  as  new  vertices,  and  then 
deleting  aU  vertices  that  lie  in  [ri,r2]. 

Label  r2  with  r. 

If  the  label  of  the  vertex  following  r2  is  (p,  c)  such  that  the  base  of  the  perpendicular  from  c 
to  Lx{t)  lies  below  r2,  then  relabel  this  vertex  with  (p,r2). 

If  ri  lies  in  an  interval  of  Tx{f)  with  label  (p',  U),  then  label  ri  with  (p',  U)  or  (p',  ri)  according 
as  the  base  of  the  perpendicular  from  U  to  Lx{t)  lies  below  or  above  vi. 

See  Figure  5. 
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Right  event:  At  a  right  event  at  time  t,  Lx{t)  touches  the  right  side,  say  [rs,  rj,  of  some 
rectangle  r  ^  R,  and  we  imagine  that  r  ceases  to  be  active  as  an  obstacle  intersecting  Lx{t).  If 
the  left  side  of  r  is  [7'i,r2],  then  r*i  =  and  r2  =  (of  course,  identifying  points  with  their 
projections  on  the  y  axis),  so  that  rs  and  r4  already  belong  to  Tx{t')  (as  they  were  inserted  at 
the  left  event  when  touched  the  left  side  of  r  and,  clearly,  could  not  have  been  removed  by 
any  intermediate  events). 

Updating  D:  Delete  rs  and  r4  from  and  assume  the  labels  of  the  intervals  just  below 

rs  and  just  above  r4  in  Lx{t')  are  (p,  c)  and  {p',c'),  respectively. 

Updating  T^'-  Three  cases  may  arise  according  as: 

(i)  r4  is  closer  to  p  then  p':  Extend  the  upper  endpoint  of  the  interval  just  below  rs  to  at 
least  r4,  and  then  proceed  to  extend  it  further  upwards  exactly  as  in  the  manner  for  a  point 
event . 

(ii)  rs  is  closer  to  p'  then  p:  Extend  the  lower  endpoint  of  the  interval  just  above  r4  to  at 
least  rs,  and  then  proceed  to  extend  it  further  downwards  exactly  as  in  the  manner  for  a  point 
event . 

(iii)  T4.  is  closer  to  and  is  closer  to  p:  Determine  the  point  r  £  (7*3,  r4)  which  is  equidistant 
from  p  and  ^ ,  and  extend  the  intervals  just  below  r3  and  just  above  r4  to  meet  at  r . 

See  Figure  5. 

Some  useful  facts  that  are  not  hard  to  verify,  and  help  bound  the  complexity  of  maintaining 
the  status  tree  Tx,  are  collected  in: 

Lemma  2  The  following  hold: 

1.  The  update  procedure  does  indeed  maintain  as  well  as  the  labels  of  its  vertices  correctly 
(i.e.,  according  to  conditions  l-f  for  the  labels  Bi  given  at  the  beginning  of  the  section). 

2.  If  the  labels  (p,  c)  and  (p',  c)  were  both  associated  with  vertices  of  Tx,  even  at  different 
times,  then  p  =  p',  implying  that  the  set  of  possible  labels  has  cardinality  0{n-\-  k). 

3.  No  two  vertices  of  Tx  at  a  given  instant  can  have  the  same  label. 

4.  If  a  label  B  is  associated  with  some  vertex  ofTx  at  time  t' ,  but  is  not  associated  with  any 
vertex  of  Tx  at  a  time  t  >  t'  (i.e.,  the  vertex  associated  with  B  has  been  deleted  by  some 
intermediate  event),  then  B  can  never  again  be  associated  with  any  vertex  of  Tx  at  a  time 
i"  >  t.  In  other  words,  once  a  vertex  with  a  given  label  is  deleted  that  label  cannot  reappear. 

5.  At  each  event,  the  update  of  Tx  involves  at  most  three  insertion  or  relabeling  operations. 
Proof. 

1 :  Follows  from  the  description  of  the  update  procedures. 

2:  If  (p,  c)  is  a  label  consider  two  cases: 

(a)  c  e  5,  in  which  case  we  must  have  p  =  c,  and, 

(b)  c  G  Cl,  in  which  case  p  is  a  closest  to  c  of  the  points  of  S  from  which  there  is  an  x-path 
to  c.  However,  after  the  choice  of  p  is  made,  at  the  time  of  updating  D  at  the  left  event  when 
Lx  lies  on  c,  it  is  never  changed. 

In  either  case,  we  see  that  (p,  c)  =  (p',  c)  p  —  p' ■ 

3:  This  follows  from  condition  4  for  the  labels  Bi  given  at  the  beginning  of  the  section. 

4 :  Follows  by  examining  each  case  where  a  deletion  may  occur  in  the  update  procedure. 

For  example,  consider  the  case  of  a  point  event  at  q  that  results  in  the  deletion  of  the  vertex, 
say  02  G  Tx  with  label  (p,  c).  This  impUes  that  the  base  c  of  the  perpendicular  from  c  to  Lx{t) 
is  closer  to  q  than  p.  See  Figure  6. 
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Figure  6:  Ulustration  for  the  proof  of  Lemma  2. 


If,  at  some  time  t"  >  t,  the  label  (p,  c)  reappears  in  T^,  say  with  vertex  a4,  then  we  must 
have  that  the  base  1  of  the  perpendicular  from  c  to  which  passes  through  c,  is  closer  to 

p  than  q.  This  is  clearly  not  possible. 

5:  Follows  by  a  scrutiny  of  each  case  of  the  update  procedure.  Q 

Lemma  2  implies: 

Lemma  3  Tx  is  of  size  0{n  +  k)  throughout,  and  the  total  number  of  operations  performed  on 
Tx  through  the  sweep  is  0{n+  k).  Therefore,  the  maintenance  of  Tx,  which  is  height-balanced, 
through  the  sweep  takes  0{{n  +  A;)log(n  +  k))  time. 

At  the  end  of  the  sweep  we  know,  from  the  contents  of  D,  for  each  point  g  €  5,  a  point  p 
which  is  closest  to  q  amongst  points  of  5,  different  from  q,  from  which  there  is  an  x-path  to  q. 
We,  therefore,  have  the  following: 

Proposition  4  In  0((n  +  A:)log(ra  +  A;))  time  one  can  determine,  following  a  sweep  of  the  plane 
in  the  x  direction,  for  each  point  q  £  S,  a  point  p  which  is  closest  to  q  amongst  points  of  S, 
different  from  q,  from  which  there  is  an  x-path  to  q. 

In  an  exactly  similar  manner  we  can  determine,  within  the  same  time  bound,  but  following 
a  sweep  of  the  plane  in  the  z  direction  (where  z  =  —x,  or  :ty),  for  each  point  q  E  S,  a  point  p 
which  is  closest  to  q  amongst  points  of  S ,  different  from  q,  from  which  there  is  a  z-path  to  q.  □ 

Now,  the  following  is  easy  to  see: 

Observation  1  For  any  two  points  p,q  €  S,  there  is  either  a  shortest  x-path  or  a  shortest 
y -path  from  one  ofp,q  to  the  other. 

Consequently,  a  closest  pair  can  be  determined  by  examining,  for  each  point  q  £  S,  points 
of  5,  different  from  q,  which  are  closest  to  q  with  distances  being  measured  to  q  either  along 
x-paths  or  y-paths.  Applying  Proposition  4  we  have: 

Theorem  1  In  0{{n  +  k)  log(n  +  k))  time  one  can  determine,  following  sweeps  of  the  plane  in 
the  X  and  y  directions,  a  closest  pair  amongst  points  of  S.  O 
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4  Finding  Nearest  Neighbors 

Determination  of  a  nearest  neighbor  for  each  point  q  6  S  \s  also  straight-forward  using  Propo¬ 
sition  4.  Relevant  is; 

Observation  2  For  any  point  p  e  S,  the  shortest  path  from  p  to  a  given  point  qe  S  is  a  z-path, 
where  z  is  either  ±x  or  ±y  (these  are,  of  course,  not  mutually  exclusive  possibilities). 

Proposition  4,  therefore,  gives: 

Theorem  2  In  0{{n  +  k)  log(n  +  k))  time  one  can  determine,  following  sweeps  of  the  plane  in 
each  of  the  four  axes-parallel  directions,  a  nearest  neighbor  for  each  point  q  E  S.  □ 

5  Computing  the  Voronoi  Diagram 
5.1  The  Plan 

Denote  the  Voronoi  diagram  of  a  point  set  S  in  the  presence  of  an  obstacle  set  R  of  rectangles 
by  Vr{S).  Due  to  the  presence  of  obstacles,  one  needs  to  be  careful  about  the  structure  of 
1/^(5):  it  is  a  planar  straight-line  graph  (PSLG)  such  that  each  face  is  either  a  rectangle  of  R, 
or  corresponds  to  one  point  of  S,  say  p,  and  is  the  locus  of  points  (£  G  —  Urgi?intertor(r))  that 
are  as  close  to  p  as  to  any  other  point  of  S.  Call  the  face  corresponding  to  p  the  V-cell  of  p,  and 
denote  it  Vr{S,p).  The  boundary  of  Vr{S,p),  denoted  Br{S,p),  is  a  polygonal  line,  consisting 
of  one  or  more  connected  components,  each  point  of  which  lies  either 

•  on  a  bisector  between  p  and  some  other  g  6  5,  or 

•  on  the  boundary  of  an  r  €  R. 

It  may  be  observed  that,  due  to  the  rectilinearity  of  the  metric,  each  straight-line  segment  of  a 
bisector  is  either  horizontal,  vertical,  or  inclined  at  45°  or  135°  to  the  positive  direction  of  the 
a:-axis.  See  Figure  7  for  illustration.  To  avoid  cumbersome  technicalities,  we  shah,  in  case  of 
bisectors  with  non-zero  area  (see  [12])  choose  vertical  lines  as  bisectors,  and,  further,  make  a 
standard  assumption  of  “general  position”  such  as  no  more  than  three  points  of  S  are  co-circular. 

The  following  bounds  the  complexity  of  Vr{S),  where  |i?|  =  n  and  |5|  =  k: 

Lemma  4  The  complexity  of  Vr(S)  as  a  planar  graph  is  0(n  -1-  k). 

Proof.  Consider  the  planar  graph  V  whose  vertex  set  consists  of 

(a)  Voronoi  vertices  (where  three  bisectors  meet),  and 

(b)  vertices  where  a  bisector  meets  a  rectangle, 
and  whose  edge  set  consists  of 

(a)  pieces  of  rectangle  boundaries,  and 

(b)  pieces  of  bisectors, 

that  join  such  vertices,  ignoring  the  individual  segments  that  comprise  each  piece. 

See  Figure  8(a)  for  the  graph  corresponding  to  the  Voronoi  diagram  of  Figure  7. 

Then  V  is  a  planar  graph  with  vertices  of  degree  three  (and  no  higher  by  the  assumption  of 
general  position),  so  by  Euler’s  formula  the  complexity  of  V  is  linear  in  the  number  of  its  faces 
which  is,  of  course,  at  most  n  +  k.  Note  that  V  may  have  parallel  edges  and  edges  going  off  to 
“infinity”.  Note  also  that  rectangles  that  do  not  intersect  bisectors  are  unrepresented  in  V  (see 
Figure  8(b)  and  (c)). 
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Figure  7:  The  Voronoi  diagram  for  a  set  S  of  seven  points  and  a  set  R  of  three  rectangles. 

The  shaded  region  is  the  V-face  Vr{S,p)  of  p.  The  separator  (dashed  line)  splits  S  into  Si  and 
S2,  and  the  bisector  between  them  is  shown  in  bold. 
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Figure  8:  (a)  The  planar  graph  V  corresponding  to  the  Voronoi  diagram  in  the  previous  figure, 
(b)  A  Voronoi  diagram  showing  the  points  of  S  UC  (indicated  in  parenthesis)  to  which  non-last 
edges  of  bisector  pieces  are  charged,  (c)  The  planar  graph  corresponding  to  the  Voronoi  diagram 


Now,  for  the  complexity  of  Vji^S)  we  have  to  take  into  account,  in  addition  to  the  complexity 
V: 

(1)  the  complexity  of  the  unrepresented  rectangles.  This  clearly  totals  0(n). 

(2)  the  complexity  of  the  pieces  of  rectangle  boundaries  that  correspond  to  edges  of  type  (a). 
Each  such  piece  comprises  at  most  4  segments,  so  the  total  complexity  is  again  0{n). 

(3)  the  complexity  of  the  pieces  of  bisectors  that  correspond  to  edges  of  type  (b).  This 
requires  more  careful  analysis  that  we  do  next. 

Analysis  of  case  (3):  We  must  count  all  segments  comprising  those  pieces  of  bisectors  that 
correspond  to  edges  of  type  (b).  For  this  we  shall  describe  an  accounting  scheme  where  ea^ 
such  segment  is  “charged”  either  to  a  point  in  S,  a  corner  of  a  rectangle  in  R,  or  an  edge  ofV. 

Consider  an  arbitrary  edge  of  V  of  type  (b).  Say  it  corresponds  to  the  piece  b  of  the  bisector 
between  points  p  and  g  of  S.  Orient  b  arbitrarily  and  imagine  proceeding  along  b,  segment  by 
segment,  from  the  first  segment  to  the  last.  Assume  each  segment  to  be  open  (i.e.,  it  does  not 
include  its  endpoints).  For  each  segment  e  of  b,  except  the  last,  denote  its  successor  by  e^ 

It  is  not  tard  to  see  that,  as  a  consequence  of  the  rectilinear  geometry,  exactly  one  of  the 
following  8  cases  must  hold,  where  each  case  is  obtained  by  replacing  m  by  p  or  g,  and  2;  by  ±a; 
or  ±y,  in  the  following  statement: 

The  shortest  path  from  m  to  each  point  of  e(e')  is  a  z-path  and  there  is  a  unique 
z-first  path  to  each  point  of  e(e'),  all  with  common  z-anchor  c(c'),  where 

(a)  either,  c', 

(b)  or,  c  =  c',  and  the  ray  in  the  z  direction  from  c  intersects  b  at  the  endpoint 
between  e  and  e'. 

(Note  that  the  z-anchor  of  a  segement  may  be  a  point  of  5  or  C,  the  set  of  corners 
of  rectangles  of  R,  and  is  unique  as  each  segment  is  open.) 

In  this  case,  charge  the  segment  e  to  the  z-anchor  c.  See  Figure  8(b). 

It  may  also  be  seen  that,  conversely,  if  a  point  in  5  U  C  does,  in  fact,  accumulate  a  charge  in 
this  case,  then  m  is  uniquely  determined  as  the  point  of  S  closest  to  c.  And,  c  can  accumulate 
at  most  2  charges  in  each  case  corresponding  to  z  =  ±x  or  ±y,  from  the  two  segments  on 
either  side  of  the  endpoint  where  the  ray  in  the  z  direction  from  c  possibly  intersects  a  bisector. 
Therefore,  each  point  of  5  U  C  may  accumulate  at  most  8  charges  under  our  accounting  scheme. 
This  proves  that  the  number  of  segments  that  are  not  the  last  segments  in  the  bisector  pieces 
to  which  they  belong  is  0(n  -f-  k).  _ 

Finally,  charge  each  last  segment  of  a  bisector  piece  to  that  edge  of  V  to  which  the  piece 
corresponds,  so  that  each  edge  of  V  is  charged  once,  proving  that  the  number  of  such  segments 
is  also  0{n  -f-  k). 

Thus,  the  number  of  aU  segments  comprising  those  pieces  of  bisectors  that  correspond  to 
edges  of  type  (b)  is  0{n  +  k),  and  this  concludes  our  analysis  of_case  (3). 

Adding  the  complexities  of  cases  (1),  (2),  and  (3)  to  that  of  V  proves  the  lemma.  □ 

The  overall  plan  for  our  Voronoi  diagram  algorithm  is  divide-and-conquer  with  an  extension 
phase: 

Apply  Proposition  3  to  find  an  a;(-i/)-path  P  such  that  the  subsets  of  S  and  R  to  the  left 
and  right  of  P  are  Si  and  Ri,  and  S2  and  R2,  respectively.  This  guarantees  that  the  number  of 
rectangles  in  Ri  and  R2  is  some  constant  fraction  of  the  number  of  rectangles  of  R  (note  that 
we  do  not  need  any  such  assumption  of  “good”  separation  on  5i  and  S2  for  our  algorithm  to 
work).  Recursively  compute  yRi(‘S'i)  and  Fr2(‘5'2)-  Then,  extend  VR.(5'j)  to  1^(5, ■),  i  =  1,2, 
and,  finally,  merge  Vr{Si)  and  Vr{S2)  with  a  Shamos-Hoey  type  scan  [16]  to  obtain  yR(5'). 
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5.2  Extending  a  Diagram 

Let  us  consider  the  problem  of  extending  Vrj(5i)  to  Vr{Si)  (extending  Vr2{S2)  to  Vr{S2)  is 
exactly  similar).  If  the  closed  region  of  the  plane  left  (right)  of  P  is  denoted  Pi  {P2),  it  is  clear 
that  Vr{Si)  r\  Pi  =  Vr^{Si)  n  Pi,  as  the  shortest  path  joining  two  points  in  Pi  lies  wholly  in 
Pi,  even  considering  all  obstacles  in  R.  Therefore,  it  remains  to  construct  V/{(5i)  0  P2,  that  is 
extend  ^^^(Pi)  right  of  P. 

Consider  points  p  6  Si,q  G  P2-  Clearly,  q  can  be  situated  in  either  the  north-west,  north¬ 
east,  or  south-east  quadrants  as  viewed  from  p,  and,  given  the  shape  of  P,  in  the  first  case  there 
can  be  no  vertical  blockage  between  them,  and  in  the  third  case  there  can  be  no  horizontal 
blockage.  Recalling  Proposition  1  we  have: 

Observation  3  A  shortest  path  from  a  point  p  E  Si  to  a  point  q  £  P2  is  either  an  x-path  or  a 
y-path. 

Therefore,  to  And  a  p  €  Pi  nearest  to  a  given  q  €  P2  (or,  equivalently,  the  V-cell  Vr{Si,p) 
containing  g),  we  need  only  consider  x-paths  and  y-paths  from  points  of  Pi  to  q,  suggesting,  in 
fact,  the  geometric  methods  of  Section  3. 

This  leads  to  the  definition  and  construction  of  z-diagrams,  for  z  =  ±a;,±y.  In  particular, 
let  us  describe  the  x-diagram. 

Constructing  the  x-diagram  of  and  R:  Sweep  the  plane  G,  containing  only  Pi  and  R, 
as  in  Section  3,  in  the  x  direction.  As  the  line  sweeps  through  the  plane,  build,  on  the  fly,  a 
planar  straight-line  graph,  called  the  i-diagram  of  Pi  with  obstacle  set  R,  and  denoted  Xr{Si), 
such  that  each  face  of  Xr{Si)  is  either 

•  a  rectangle  of  R,  or 

•  corresponds  to  one  point  of  Pi  X  (Pi  U  Ci),  say  (p,c)  (where  Ci  is  the  set  of  corners  of 
left  edges  of  rectangles  of  P;  see  the  discussion  of  labels  at  the  beginning  of  Section  3),  in 
which  case  it  is  the  locus  of  those  points  g  such  that,  amongst  the  points  of  Pi  from  which 
there  is  an  x-path  to  g,  p  is  the  nearest  and  c  is  the  x-anchor  of  g  w.r.t.  p,  or 

•  the  one  remaining  face,  not  of  the  preceding  two  types,  that  is  the  locus  of  those  points 
q£  G  -  Dr^R  interior{r)  that  cannot  be  reached  from  any  point  of  Pi  by  an  x-path. 

See  Figure  9  for  illustration.  Xr{Si)  can  be  built,  and  represented  as  a  doubly-connected-edge- 
list  (DCEL,  see  [15]),  by  tracing  the  motion  on  the  plane  of  points  of  the  status  tree  (lying 
on  Lxj  of  course:  imagine  an  inkspot  at  each  such  point),  and,  when  the  sweep-line  Lx  touches  a 
point  p  G  Pi  (or,  c  G  Ci)  at  time  t,  tracing  the  motion  on  the  plane  (along  Lx{t))  as  we  proceed 
in  either  direction  from  p  (or,  possibly,  c)  to  determine  the  extent  of  the  new  interval  with  label 
{p,p)  (or,  possibly,  (p,c):  imagine  marking  this  new  interval  with  a  pen;  see  the  discussion  on 
updating  the  status  tree  in  Section  3). 

Observation  4  Each  face  of  Xr{Si)  is  a  rectilinear  polygon,  and  the  face  of  Xr{Si)  corre¬ 
sponding  to  (p,  c)  G  Pi  X  (Pi  U  Cl),  called  the  x-ceU  of  (p,  c)  and  denoted  Xr{Si,  (p,  c)),  contains 
c  on  its  boundary.  Xr{Si)  is,  of  course,  a  planar  graph  analogous  to  a  Voronoi  diagram  Vr{Si), 
but  where  distances  are  measured  only  along  x-paths  and  that,  further,  has  been  refined  up  to 
anchors. 

Assume  the  number  of  points  in  P,-,  i  =  1,2,  is  k;.  The  following  bounds  the  complexity  of 
the  x-diagram,  as  well  as  that  of  computing  it: 
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Figure  9;  The  x-diagram  Xr{Si),  for  5i  =  {pi,P2,P3,P4,P5,P6,P7}  and  R  = 

{7-1,72,73,74,75,76,77}.  The  shaded  region  m  the  left  consists  of  points  that  cannot  be 
reached  from  points  of  Si  by  x-paths,  while  tne  shaded  region  in  the  center,  for  example,  is 


Lemma  5  The  complexity  of  X{i{S-i)  as  a  planar  graph  is  0{n+  ki),  and  the  time  to  compute 
it  is  0{{n  +  A;i)log(n  +  ki)). 

Proof.  The  complexity  of  XftfSi)  as  a  planar  graph  follows  from  standard  geometric  arguments, 
while  the  time  to  compute  it  is  dominated  by  the  time  to  sweep  the  plane  with  L^,  which  is 
0((n  +  fci)log(n  + /ji))  from  Lemma  3.  ° 

This  completes  a  description  of  constructing  the  a:-diagram 

One  can  define  and  construct,  in  an  exactly  similar  manner,  the  z-diagrams  Zfi(Si),  for 
z  =  -x,±y  {Z  =  -X,±Y,  respectively).  In  particular,  construct  the  y-diagram  Yr{Si). 

Now  note  that,  by  Observation  3,  Xr(S\)  0  P2  and  0  P2  together  contain  essentially 

aU  the  proximity  information  of  VR{Si)r\P2.  How  then  do  we  use  these  two  diagram  to  explicitly 
construct  Vr{Si)  D  P2?  Our  plan  is  to  scan  the  plane  upwards  with  a  horizontal  sweep-line  Ly, 
starting  at  a  position  sufficiently  far  south  and  constructing  the  part  of  yR(5i)  0  P2  beneath  Ly 
as  we  proceed. 

Before  describing  this  construction,  we  record  certain  relevant  geometric  features  of  VH(5'i)n 

P2- 

Proposition  5  The  following  hold: 

1.  Each  edge  e  ofVR{Si)  fl  P2  is  a  straight-line  segment  that  is  either  horizontal,  vertical,  or 
inclined  ai  45°  to  the  positive  direction  of  the  x-axis.  In  particular,  no  edge  ofVR(Si)nP2 
is  inclined  at  135°  to  the  positive  direction  of  the  x-axis. 

2.  Each  edge  e  ofVn{Si)  0  P2  is  one  of  the  following  five  mutually  exclusive  types: 

•  P-edge:  e  is  part  of  an  edge  of  P. 

•  R-edge:  e  is  part  of  an  edge  of  a  rectangle  of  R  within  P2- 

•  XY-edge:  e  is  part  of  the  bisector  between  two  points  p,q  G  Si,  such  that  the  shortest 
path  to  any  point  of  e  from  p  is  an  x-path  and  the  shortest  path  to  any  point  of  e  from 
q  is  a  y-path. 

•  XX -edge:  e  is  part  of  the  bisector  between  two  points  p,q  £  Si,  such  that  e  is  not  of 
type  XY  and  the  shortest  paths  to  any  point  of  e  from  both  p  and  q  are  x -paths. 

•  YY -edge:  e  is  part  of  the  bisector  between  two  points  p,q  £  Si,  such  that  e  is  not  of 
type  XY  and  the  shortest  paths  to  any  point  of  e  from  both  p  and  q  are  y-paths. 

3.  An  inclined  edge  is  always  an  XY-edge,  a  vertical  edge  may  be  either  an  XY-  or  YY-edge, 
while  a  horizontal  edge  may  be  either  an  XY-  or  XX -edge. 

4.  Each  vertex  v  0/ Vr(5i)  H  P2  is  one  of  the  following  six  mutually  exclusive  types: 

•  P-vertex:  v  is  a  corner  of  P. 

•  R-vertex:  v  is  a  corner  of  a  rectangle  of  R  within  P2. 

•  B -vertex:  v  is  a  point  at  which  the  bisector  between  two  points  of  Si  turns  (from  one 
straight-line  segment  to  another)  within  P2* 

•  BP -vertex:  v  is  the  intersection  with  P  of  the  bisector  between  two  points  of  Si. 

•  B  R-vertex:  v  is  the  intersection  with  an  edge  of  a  rectangle  of  R,  within  P2,  of  the 
bisector  between  two  points  of  Si . 
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•  V -vertex:  v  is  the  intersection  of  three  bisectors  within  (i.e.,  a  Voronoi  vertex). 

Proof.  The  classification  of  the  edges  and  vertices  of  Vii{Si)  fl  P2  (items  2  and  4)  follows 
immediately  from  the  definition  of  Vr{S\)  and  Observation  3.  (See  Figure  10  for  illustration.) 

For  the  proof  of  item  1,  observe  that,  given  two  points  pi  and  p2  such  that  a  segment  of  the 
bisector  between  them  is  inclined  at  135°  to  the  positive  direction  of  the  x-axis,  and  given  a 
subsegment  d  of  that  inclined  segment,  one  of  pi  and  p2  must  lie  in  the  upper  right  “quadrant” 
determined  by  d  (refer  Figure  11(a)).  Now,  suppose  that  there  exists  an  edge  e  of  Vr(S'i)  n  P2 
inclined  at  135°  to  the  positive  direction  of  the  x-axis.  This  would  then  imply  (refer  Figure  11(b)) 
that  a  point  of  S'!  lies  to  the  right  of  P,  which  is  impossible. 

For  item  3,  we  can  use  an  argument  similar  to  the  one  given  above  to  show  that  an  XX-edge 
can  run  only  horizontally,  and  a  YY -edge  only  vertically.  We  leave  the  details  to  the  reader.  □ 

Observe  that  at  any  instant  through  the  intended  sweep,  when  Ly  lies  along  Ly{t),  Ly{t)  0 
(V^(5i)  n  P2)  consists  of  zero  or  more  horizontal  edges  of  Vr(Si)  H  P2,  and  the  intersections  of 
some  vertical  and  inclined  edges  of  Vr(5i)  (1  P2  with  Ly{t)  (these  intersections  may,  in  fact,  be 
vertices  of  Vii(5i)  n  P2).  We  keep  a  description  of  Ly(t)  D  (Vft(5i)  H  P2)  as  a  sorted  sequence 
{xj  :  i  =  0, . . .}  of  points  stored  in  records  of  a  height-balanced  status  tree  Ty.  Each  x,-  is  either 
the  endpoint  of  a  horizontal  edge  of  yR(5i)  n  P2,  or  the  intersection  of  a  vertical  or  inclined 
edge  of  ^^(^i)  n  P2  with  Ly{t),  which  information  is  stored  at  the  record  for  x,-  (of  course, 
XO  =  P  n  Ly{t)). 

Our  plan  is  to  mcdntain  Ty  (at  least  impbcitly,  but  not  necessarily  exactly,  which  we  justify 
later)  through  the  sweep.  This  will  allow  us  to  build  V’k(5i)  n  P2,  and  represent  it  as  a  DCEL, 
on  the  fly. 

Initialize  the  event  point  schedule  E  by  including  in  it  the  following  points  sorted  by  y- 
coordinate: 

1.  P- vertices  of  V/?(5i)  0  P2  that  are  left  endpoints  of  horizontal  edges  of  P. 

2.  P- vertices  of  yR(5'i)  H  P2  that  are  left  endpoints  of  horizontal  edges  of  r  G  ^2- 

3.  PP- vertices  of  Vfl(5i)  n  P2. 

4.  PP- vertices  of  Vr(5i)  n  P2  that  lie  on  upper  or  right  edges  of  rectangles  of  P2. 

Points  of  the  third  and  fourth  type  are  called  start  vertices  as  they  mark  “starting”  endpoints 
for  bisectors  in  Vr{Si)  H  P2  . 

It  should  be  remarked  that  new  event  points,  either  so-called  PP-vertices  (to  be  defined)  or 
F- vertices,  will  be  added  to  E  as  the  sweep  proceeds.  F-vertices  are  also  start  vertices. 

We  claim: 

Lemma  6  E  can  be  initialized  in  0{{n-\-  A:i)log(n-|-  kf))  time,  and  within  that  same  time  one 
can  determine,  for  each  start  vertex  v  of  E,  the  direction  and  type  of  the  edge,  say  e,  that  starts 

at  that  vertex,  as  well  as  the  identity  of  the  bisector  of  which  e  is  a  part  (i.e.,  the  two  points  of 

Si  that  this  bisector,  in  fact,  bisects). 

Proof  It  is  trivial  to  include  points  of  the  first  two  types.  Points  of  the  third  type  can  be 
determined  as  the  intersections  of  edges  of  Vkj(5i)  (which  has  been  recursively  computed)  with 

P. 

For  points  of  the  fourth  type,  we  claim  that,  in  fact,  aU  PP-vertices  of  Vr{Si)  n  P2  can  be 
determined  in  0((n  -t-  A:i)log(n  4-  ki))  time  by  the  following  steps: 
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Figure  11;  Illustration  for  the  proof  of  Proposition  5:  (a)  The  disposition  of  two  points  whose 
bisector  has  a  segment  inclined  at  135°  to  the  positive  direction  of  the  i-axis,  and  (b)  part  of 
such  a  bisector  to  the  right  of  P. 

(a)  Preprocess  the  planar  graphs  XhC*?!)  and  yR(5i),  for  logarithmic  time  point  location,  in 
time  0((n  +  /:i)log(n  +  ^i))  (see  [15]). 

(b)  Locate  the  lower  left  corner  of  each  rectangle  of  R2  in  both  graphs  Xr(Si)  and  yH(5i) 
in  total  time  0{nlog{n  +  ki)). 

(c)  Traverse  clockwise  the  boundary  of  each  r  G  R2,  starting  from  its  lower  left  corner,  while 
tracking,  at  the  same  time,  the  nearest  points  of  and  S2-  For  this,  by  Observation  3,  we  sim¬ 
ply  need  to  determine  the  successive  x-cells  and  j/-cells  that  are  crossed  through  the  traversal. 
And,  during  the  traversal,  a  point  v  on  the  boundary  of  r  that  is  nearest  to  two  distinct  points 
of  Si  is  marked  as  a  5  i?- vertex. 

Before  describing  the  rather  complex  procedure  for  updating  Ty  at  each  event  point,  let  us 
try  to  provide  the  guiding  geometric  insight  which  is  relatively  simple. 

View  n  F2  as  a  “set  of  (connected)  bisector- components”,  each  bisector-component 

starting  at  either  a  5P- vertex  or  a  5 P- vertex  (on  the  upper  or  right  edge  of  a  rectangle  of 
R),  and  traveling  monotonicaUy  in  the  upwards  direction  either  to  “infinity”,  or  ending  at  a 
V-vertex  or  on  the  lower  or  left  edge  of  a  rectangle.  This  view  is  motivated  by  the  first  item  of 
Proposition  5,  and  justifies  the  inclusion  of  the  third  and  fourth  types  of  points  in  the  schedule 
E  (of  course,  points  of  the  first  and  second  types  are  included  in  E  as  they  mark  the  predictable 
events  of  including  edges  of  P  and  rectangles  of  R  into  Vr(Si)  D  P2). 

Now,  a  yy-bisector-component  (i.e.,  a  bisector-component  consisting  of  yy-edges:  it  is  not 
hard  to  see  that  if  one  edge  of  a  bisector  component  is  of  type  yy  {XY,  XX),  then  every  edge 
of  that  bisector-component  is  of  type  yy  (XY,  XX))  travels  vertically  upwards  from  its  start 
vertex,  either  indefinitely,  or  till  it  ends  either  at  a  V-vertex,  or  on  the  lower  edge  of  a  rectangle 
(the  second  case  may  be  viewed  as  the  yy-bisector-component  being  struck  from  the  left  by  a 
horizontal  edge  of  Vr(5i)  H  P2)- 

An  XX-bisector-component  travels  horizontally  rightwards  from  its  start  vertex,  either  in¬ 
definitely,  or  tiU  it  ends  at  a  V-vertex,  or  on  the  left  edge  of  a  rectangle  (the  second  case  may 
be  viewed  as  the  XX-bisector-component  striking  a  vertical  edge  of  Vr(5i)  fl  P2). 

With  an  Xy-bisector-component  we  must  be  more  careful  as  it  may  make  turns  while  trav¬ 
eling  from  its  start  vertex.  Even  so,  the  Xy-bisector-component  may  be  tracked  with  the  help 
of  the  X-  and  y-diagrams,  for  within  a  fixed  a;-cell  and  a  fixed  y-cell  (where  the  x-anchors  and 
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y-anchors  of  points  w.r.t.  nearest  points  in  do  not  change)  the  bisector-component  makes 
at  most  two  turns  and  its  path  may  be  determined  easily.  Thus,  we  track  the  Xy-bisector- 
component  first  from  its  start  vertex  n,  through  the  intersection  of  the  x-ceU  and  r/-ceU  to  which 
V  belongs,  to  the  point,  say  where  it  first  leaves  one  of  these  two  cells  to  enter  a  new  x- 
or  y-cell.  Note  that  the  Xy-bisector-component  must  make  a  turn  at  vi,  and  we  repeat  the 
tracking  procedure  starting  from  vi,  . . .,  either  indefinitely,  or  till  it  ends  either  at  a  V-vertex 
or  on  the  lower  or  left  edge  of  a  rectangle. 

With  this  insight  in  mind,  we  shall  describe  the  detailed  procedure  for  updating  Ty.  However, 
before  beginning  the  procedure  we  need  a  preprocessing  stage: 

Preprocess  each  x-ceU  of  X/e(5i)  and  y-ceU  of  for  O(log(n  -1-  A;i))-time  ray-shooting 

(i.e.,  to  determine  where  a  ray  fired  in  some  given  direction  from  some  given  point  in  that  cell 
wiU  first  strike  the  boundary  of  the  cell)  in  total  0{n  4-  ki)  time  (see  [9]  for  linear-time  prepro¬ 
cessing  for  ray-shooting).  Further,  preprocess  the  horizontal  edges  of  both  P  and  rectangles  of 
R,  again  in  total  0{n  +  ki)  time,  for  0(log(n  +  fci))-time  vertical  ray-shooting  (i.e.,  determining 
which  of  these  horizontal  edges,  if  any,  will  be  struck  by  a  vertical  ray  fired  from  some  given 
point  on  the  plane,  see  the  computation  of  vertical  adjacency  maps  [15,  p.  349]). 

Procedure  for  updating  Ty:  Initially,  when  Ly{t)  lies  south  of  any  point  in  the  schedule 
E,  Ty{t)  contains  only  one  point,  the  intersection  of  Ly{t)  with  the  lowest  vertical  edge  of  P . 
Next,  start  scanning  upwards  with  Ly,  stopping  at  successive  event  points  of  E  to  update  Ty  as 
follows: 

1.  At  a  P-vertex  v,  v  becomes  the  new  left  endpoint  of  the  intersection  Ly{t)  n  (yR(*S'i)  H  J’z)- 
Insert  v  “appropriately”  into  Ty{t)  (this  includes  recording  the  horizontal  and  vertical 
edges  of  P  in  Vr(5i)  n  P2  that  start  at  v),  and  delete  the  previous  leftmost  point,  say  ui, 
of  Ty  (this  represents  that  a  vertical  edge  of  Vr(5i)  n  P2  ends  at  ui). 

2.  At  an  P-vertex  that  is  on  the  lower  edge  of  an  r  e  P,  appropriately  insert  both  endpoints 
of  that  lower  edge  into  Ty{t)  (the  lower  edge  itself,  as  weU  as  the  vertical  edges  starting  at 
either  endpoint  of  that  edge  are  parts  of  V’r(5i)  D  P2). 

At  an  P-vertex  that  is  on  the  upper  edge  of  an  r  6  P,  appropriately  delete  both  endpoints 
of  that  edge  from  Ty{t)  (they  both  mark  endpoints  of  verticals  edges  of  yR(5i)  n  P2,  while 
the  horizontal  edge  joining  them  is  part  of  Vr(5i)  n  P2). 

3.  At  a  start  vertex  v  (i.e.,  either  a  BP-,  BR-,  or  V-vertex),  insert  it  into  Ty{t)  as  the  starting 
endpoint  of  an  edge  e  of  Vr(5'i)  n  P2,  and  also  record  the  direction  and  type  of  e,  as  weU 
as  the  identity  of  the  bisector  of  which  e  is  a  part  (see  Lemma  6). 

If  e  is  a  yy-edge,  record  the  point,  say  vi  (possibly  -foo),  at  which  the  ray  along  e  (i.e., 
vertically  upwards)  first  hits  a  rectangle  of  P,  as  the  “likely  endpoint”  of  e  (preprocessing 
for  vertical  ray-shooting  aUows  us  to  compute  vi  in  time  0(log(n  -|-  ki))),  by  inserting  v\ 
appropriately  into  E  as  an  LE-vertex. 

If  e  is  an  XX-edge,  it  is  horizontal  and,  in  fact,  lies  whoUy  along  Ly{t)  and  ends  at  the 
next  vertex  (to  the  right),  say  vi,  on  Ty{t)  which  must  correspond  either  (a)  to  a  vertical 
edge  lying  on  the  left  side  a  rectangle  of  r  G  P,  or  (b)  a  vertical  edge  that  is  part  of 
some  bisector.  In  both  cases  the  vertical  edge  ends  at  vi  (so,  also  remove  the  TP-vertex 
corresponding  to  this  vertical  edge  from  E).  In  case  (a),  a  “new”  vertical  edge  continues 
up  the  left  side  of  r  starting  from  v\.  In  case  (b),  insert  vi  into  E  a-s  a,  new  V-vertex 
which  is  the  start  of  an  Xy-edge  e'  that  is  part  of  a  bisector  whose  identity  is  also  triviaUy 
determined. 
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If  e  is  an  XY-edge,  extend  the  bisector-component  6  of  which  it  is  a  part  to  the  point, 
say  ni,  where  b  first  strikes  the  boundary  of  either  the  x-ceU  or  y-ceU  to  which  v  belongs 
(determining  vi  takes  time  0(log(n  +  ^i))  as  we  need  to  foUow  the  bisector  through  at 
most  two  turns  and,  possibly,  perform  ray-shooting).  Insert  vi  into  E  as  an  L£^-vertex, 
recording  it  as  the  ‘^likely  endpoint”  of  b, 

4.  At  an  vertex  vi  that  lies  on  the  edge  (either  left  or  lower)  of  a  rectangle,  mark  vi  as 
the  actual  endpoint  of  the  corresponding  edge  of  a  bisector-component  of  VH(5i)n  P2-  At 
an  vertex  vi  where  an  XY-edge  strikes  the  boundary  of  either  an  a;-ceU  or  a  j/'CeU, 
again  mark  vi  as  the  endpoint  of  the  corresponding  edge  of  a  bisector-component  b  of 
Vr(*S'i)  n  P2*  It  is,  in  fact,  a  turn  point  of  b.  Therefore,  as  in  the  last  paragraph  of  Step  3, 
further  extend  b  to  the  point,  say  V2,  where  it  next  strikes  the  boundary  of  either  the  a;-ceU 
or  y-cell  to  which  vi  belongs  (again,  determining  V2  takes  time  0(log(n  -f  ki))).  Insert  V2 
into  E  as  an  TjE-vertex,  recording  it  as  the  “likely  endpoint”  of 

Comment  The  points  in  Ty{t)  may  not  exactly  represent  the  intersections,  at  each  instant  t,  of 
Ly{t)  with  yR(5i)  n  P2j  because  of  the  inclined  edges  of  Vr{Si)  H  P2-  However,  this  does  not 
compromise  either  the  procedure  for  updating  Ty  or  constructing  the  DCEL  for  Vr{Si)  PI  P2, 
as  through  the  sweep  we  do  succeed  in  locating  both  endpoints  of  every  edge  of  Vr{Si)  D  P2, 

i.e.,  we  do  not  miss  any  intersections  of  edges  even  though  they  may  be  traveling  at  an  incline 
between  event  points  (in  fact,  each  intersection  wiU  be  discovered  as  an  event  point). 

It  may  be  checked  that  there  will  be  at  most  0{n  +  ki)  updates  to  Ty  at  a  cost  of  O(log(n  + 
ki))  per  update  (the  easiest  way  to  count  updates  is  to  see  the  correspondence  between  each 
update  and  a  vertex  of  Vk(5i)  fl  P2).  The  only  relevant  non-trivial  observation  here  is  that  each 
O(log(n  +  ki))  cost  to  track  an  Xy-bisector-component  b  from  one  PP- vertex  v  to  another  may 
be  “charged”  to  v  (note  that  b  must  turn  at  v,  so  that  v  is,  indeed,  a  vertex  of  Vr{Si)  H  P2), 
implying  that  the  total  cost  of  tracking  aU  Xy-bisector-components  is  0{{n+ki)log{n-{-ki)),  (A 
simple  intuition,  in  fact,  a  motivation  for  the  entire  procedure,  is  that,  to  trace  an  XX-  or  yy- 
bisector-component  we  “pay”  only  a  constant  amount,  while,  for  an  Xy-bisector-component, 
we  “pay”  only  when  it  turns.)  We  have,  therefore: 

Proposition  6  We  can  extend  VR-{Si)  to  VR^Si),  i  =  1,2,  in  0((n  +  k)\og{n+  k))  time,  □ 

5.3  Merging  Two  Diagrams 

Now,  consider  the  problem  of  finding  Vr(5),  given  Vr{Si)  and  Vr{S2)-  We  intend  a  Shamos- 
Hoey  type  scan  ([16],  see  also  [15])  to  compute  the  bisector  of  Si  and  82^  denoted  6(5'i,5'2)> 
consisting  of  those  points  q  G  —  interioT{r)),  such  that  if  pi  and  p2  are  the  nearest  to 
q  amongst  points  of  Si  and  52,  resp.,  then  d{q,pi)  =  d{q,p2).  But  first: 

Lemma  7  The  following  hold: 

1.  The  bisector  b{Si,S2)  is  a  subgraph  ofVR{S)  and,  therefore,  has  complexity  0(n  +  k), 

2.  It  consists  of  components,  each  a  polygonal  line  that,  in  either  direction,  is  either  unbounded 
or  ends  on  the  boundary  of  an  r  ^  R. 

3.  It  divides  G  —  Ur^Rr  into  two  (not  necessarily  connected)  subsets  6“(5i,  52)  and  &‘^(5i,  52) 
to  its  left  and  right,  resp,,  w.r.t,  to  some  total  orientation, 

I 

VR{S)  =  {VRiS^)nb-{Si,S2))  U  KSi,S2)  [j  {VR{S2)nb+{Si,S2))  U  ^renr. 
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Proof.  Follows  by  standard  geometric  arguments,  but  see  Figure  7.  O 

For  the  computation  of  b{Si,S2)  it  is  observed  that: 

The  bisector  b{Si,  S2)  may  be  found  by  a  scan  exactly  similar  to  the  Shamos-Hoey  type  scan 
for  an  Euclidean  space  without  obstacles,  except 

(a)  that  b{Si,S2)  may  have  turns  even  within  one  given  F-ceU  because  of  the  nature  of  L\- 
bisectors,  and  because  it  passes  through  regions  with  different  anchors  w.r.t.  the  closest  points 
in  Si  or  S2, 

(b)  that  such  turns  may  be  found  by  simultaneously  tracking  b{Si,S2)  through  the  2^- 
diagrams  ZR{Si)  {Z  =  ±X,±Y,  i  =  1,2),  as  the  2r-cells  to  which  a  point  p  G  b{Si,S2)  belongs 
determine  its  anchors  w.r.t.  the  closest  points  in  and  S2,  and,  therefore,  allow  0(log(n  +  k)) 
time  determination  of  the  straight  segment  of  6(5i,  52)  on  which  p  lies;  of  course,  this  requires 
0{{n  +  A;)log(n  +  k))  time  preprocessing  of  the  z-diagrams,  see  the  analogous  procedure  for 
tracking  an  XT-bisector-component  in  the  procedure  for  updating  Ty  in  the  previous  section, 
and, 

(c)  that  the  number  of  such  turns  is  0(n  +  fc). 

We  have,  therefore: 

Proposition  7  We  can  merge  Vi?(5'i)  and  Vr{S2)  to  obtain  Vr{S)  in  0{{n  +  k)log{n  +  A;)) 
time. 

5.4  Putting  Everything  Together 

We  now  have  aU  the  pieces  required  to  execute  the  plan,  described  in  Section  5.1,  for  an  algorithm 
computing  the  Voronoi  diagram: 

Given  the  problem  of  computing  the  Voronoi  diagram  defined  by  S  and  R,  where  |51  =  fc  and 
|i2|  =  n,  we  first  subdivide  (using  Proposition  3)  into  two  subproblems  defined  by  5i  and  i?i, 
and  S2  and  R2,  resp.,  where  1511  =  ki,  |i?i|  =  ni,  |52|  =  k2,  |i?2|  =  ^2,  and  7z.i,n2  >  |n.  This 
subdivision  requires  0((n+A:)logn)  time,  as  it  takes  0(nlog  n)  time  to  build  the  the  subdividing 
path  P  and,  subsequently,  0(log7i)  time  per  element  of  5  U  jR  to  determine  on  which  side  of  P 
it  lies  by  a  binary  search  over  the  0(n)-sized  P. 

After  recursively  finding  the  diagrams  that  solve  these  two  subproblems,  we  extend  and  then 
merge  these  two  diagrams  (using  Propositions  6  and  7,  resp.)  to  obtain  a  solution  to  the  original 
problem.  These  two  phases  require  a  total  time  of  0((n  +  A;)log(n  +  k)). 

Observing  that  «dl  the  base  cases  of  the  recursion  can  be  solved  in  a  total  time  of  0{k\ogk) 
(each  base  case  is,  in  fact,  defined  by  some  number,  say  ki,  of  points  of  S  and  zero  or  one 
rectangles  of  R,  and  can  be  solved  in  0{ki\ogki)  time  by  a  modification  of  the  algorithm  of 
[12]),  it  may  be  checked  that  a  recurrence  bounds  the  running  time  of  the  entire  algorithm  to 
0{{n  +  fc)log(n  +  A;)logn). 

We  have,  therefore: 

Theorem  3  In  0{{n  + k)log{n  +  k)  log  n)  time  one  can  determine  the  Voronoi  diagram  of  5.0 

6  Computing  a  Minimal  Spanning  Tree 

The  problem  of  computing  a  minimal  spanning  tree  of  S  (in  the  presence  of  the  obstacle  set  R) 
may,  in  fact,  be  reduced  in  linear  time  to  that  of  computing  the  Voronoi  diagram  Vh(5).  The 
relevant  observation  is  (see  [15]  for  the  insight): 
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Observation  5  Computing  a  minimal  spanning  tree  of  S  is  equivalent  to  computing  the  minimal 
spanning  tree  of  the  dual  graph  of  Vr{S),  where  the  weight  of  an  edge  joining  the  two  vertices 
that  are  dual  to  the  faces  Vfi{S,p)  and  Vh(5,  g),  resp.,  ofVpfS),  is  d{p,q). 

Since  the  dual  graph  of  Vr{S)  is  a  planar  graph  of  size  0{n  +  k)  (and  may  also  be  determined 
from  Vr{S)  in  time  0{n  +  k)),  the  algorithm  of  Cheriton-Tarjan  [4]  finds  a  minimal  spanning 
tree  in  0(n  +  k)  time.  This  completes  the  linear  time  reduction  of  the  minimal  spanning  tree 
problem  to  the  Voronoi  diagram  problem,  and  gives; 

Theorem  4  In  0((n  +  k)  log(n  +  k)  log  n)  time  one  can  determine  a  minimal  spanning  tree  of 
S.  ° 

7  Conclusions 

For  each  of  the  four  problems  (a)-(d)  described  in  the  introduction,  reduction  from  the  element 
uniqueness  problem  |7]  provides  an  Cl^klo^k)  time  lower  bound,  while  an  Q,(^n  -t-  k)  time  lower 
bound  holds  trivially  because  of  the  input  size.  This  gives  the  f2(n  +  klogk)  time  lower  bound 
which,  to  our  knowledge,  is  the  best  available,  and  shows  that  each  of  our  four  algorithms  is 
near-optimal.  Whether  they  can  be  improved  is,  of  course,  open. 

However,  it  seems  unlikely  that  the  time  wOl  be  lowered  below  C)((n-|-A:)  log(n-l-A:))  as,  in  some 
sense,  each  of  the  four  problems  has  an  intrinsic  complexity  of  Cl{{n  +  k)\og{n  +  k)):  for,  if  a  part 
of  each  problem  were  to  determine  the  validity  of  the  input  sets  (i.e.,  aU  rectangles  are  disjoint 
and  all  points  lie  outside  the  rectangles)  then,  indeed,  a  reduction  from  element  uniqueness  (by 
considering  point-sized  rectangles)  implies  an  fl((n  fc)log(n  +  k))  time  complexity. 

We  mention  that  one  apparently  feasible  approach  to  obtaining  0{{n  4-  A;)log(n  -f  k))  time 
complexity  -  avoid  divide-and-conquer  by  first  computing  all  four  z-diagrams  Zr{S),  Z  - 
±X,  ±y,  and  then  construct  the  entire  Voronoi  diagram  Vr{S)  with  one  plane-sweep  similar 
to  that  in  Section  5.2  -  runs  into  difficulties  that  we  believe  cannot  be  overcome  in  0((n  -h 
jb)log(n  -1-  k))  time.  For,  recall  that,  to  successfully  complete  the  sweep  of  Section  5.2  in 
0((n  +  A:i)log(n  +  ki))  time,  it  was  crucial  that  an  XX-  {YY-,  resp.)  bisector-component 
simply  runs  horizontally  (vertically,  resp.)  until  either  it  hits  a  rectangle  or  is  hit  from  below 
(left,  resp.)  by  another  Voronoi  edge.  Hence,  processing  (at  0(1)  cost)  for  a  bisector-component 
when  it  enters  a  new  x-  or  j/-cell  was  necessary  only  if  it  was  an  .W-bisector-component  (and 
an  XY- bisector-component  does,  in  fact,  make  a  turn  at  every  such  entrance,  so  “justifying” 
the  0(1)  cost).  However,  it  seems  an  analogous  property  does  not  hold  if  we  try  to  construct 
the  entire  Voronoi  diagram  in  a  single  sweep  using  the  four  z-diagrams,  and  it  is  not  clear  how 
to  avoid  expending  0(1)  processing  cost  each  time  any  bisector-component  enters  a  new  x-  or 
j/-ceU. 

Other  directions  to  consider  include  extending  or  modifying  the  techniques  described  here 
to  deal  with  more  general  classes  of  obstacles  and  different  metrics.  For  example,  we  believe 
that  similar  techniques  will  work  with  convex  polygonal  obstacles,  as  well  as  “fixed  orientation 
metrics”  [17]. 
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